;=====================================================
load "myfun_dotshade_con_vec_pattern.ncl"
load "myfun_shade_con_vector_pattern.ncl"
;=====================================================
begin

   idir = (/ "ENSO_20210719/pnmean/", "ENSO_20210719/", "ENSO_20210719/ttest/ems/" /)
   fileNames = (/"EM_rmap_ALL_atn.cdf", "EM_rmap_ALL_nino34.cdf", \
                 "EM_diff_ALL_atn.cdf", "EM_diff_ALL_nino34.cdf", \
    "HadSST_corr_n_rmap_atn_ts_4seasons.cdf", "HadSST_corr_n_rmap_nino34_ts_4seasons.cdf", \
    "HadSST_corr_n_rmap_ATN_UAS_ttest.nc", "HadSST_corr_n_rmap_NINO34_UAS_ttest.nc", \
    "HadSST_corr_n_rmap_ATN_VAS_ttest.nc", "HadSST_corr_n_rmap_NINO34_VAS_ttest.nc" /)

;=== create mask ===;
   mskfile = addfile("ENSO_20210719/remap/HadSST_corr_n_rmap_atn_ts_4seasons.remap.nc","r")
   lndmsk = mskfile->RMAP_ATNJJA0_V_TSD0JF1
   lndmsk = lndmsk - lndmsk
   tsmk1file = addfile("ENSO_20210719/ttest/ems/HadSST_corr_n_rmap_ATN_TS_ttest.nc","r")
   tsmk1     = tsmk1file->PVAL_RMAP_ATNJJA0_V_TSD0JF1
   tsmk1     = where(tsmk1.ge.0, 5, -5)
   tsmk2file = addfile("ENSO_20210719/ttest/ems/HadSST_corr_n_rmap_NINO34_TS_ttest.nc","r")
   tsmk2     = tsmk2file->PVAL_RMAP_NINO34JJA0_V_TSJJA0
   tsmk2     = where(tsmk2.ge.0, 5, -5)
;=== first data ===;
   i1file = addfile(idir(0)+fileNames(0),"r")
   i2file = addfile(idir(0)+fileNames(1),"r")
   im1file = addfile(idir(0)+fileNames(2),"r")
   im2file = addfile(idir(0)+fileNames(3),"r")
   sl1 = i1file->BX944EM_RMAP_ATNJJA0_V_TSD0JF1
   ul1 = i1file->BX944EM_RMAP_ATNJJA0_V_UASD0JF1
   vl1 = i1file->BX944EM_RMAP_ATNJJA0_V_VASD0JF1
   sr1 = i2file->BX944EM_RMAP_NINO34JJA0_V_TSJJA0
   ur1 = i2file->BX944EM_RMAP_NINO34JJA0_V_UASJJA0
   vr1 = i2file->BX944EM_RMAP_NINO34JJA0_V_VASJJA0
   sml1 = im1file->BX944EM_DIFF_ATNJJA0_V_TSD0JF1
   uml1 = im1file->BX944EM_DIFF_ATNJJA0_V_UASD0JF1
   vml1 = im1file->BX944EM_DIFF_ATNJJA0_V_VASD0JF1
   smr1 = im2file->BX944EM_DIFF_NINO34JJA0_V_TSJJA0
   umr1 = im2file->BX944EM_DIFF_NINO34JJA0_V_UASJJA0
   vmr1 = im2file->BX944EM_DIFF_NINO34JJA0_V_VASJJA0
   sml1 = sml1 + lndmsk
   uml1 = uml1 + lndmsk
   vml1 = vml1 + lndmsk
   smr1 = smr1 + lndmsk
   umr1 = umr1 + lndmsk
   vmr1 = vmr1 + lndmsk
;=== second data ===;
   sl2 = i1file->BY473EM_RMAP_ATNJJA0_V_TSD0JF1
   ul2 = i1file->BY473EM_RMAP_ATNJJA0_V_UASD0JF1
   vl2 = i1file->BY473EM_RMAP_ATNJJA0_V_VASD0JF1
   sr2 = i2file->BY473EM_RMAP_NINO34JJA0_V_TSJJA0
   ur2 = i2file->BY473EM_RMAP_NINO34JJA0_V_UASJJA0
   vr2 = i2file->BY473EM_RMAP_NINO34JJA0_V_VASJJA0
   sml2 = im1file->BY473EM_DIFF_ATNJJA0_V_TSD0JF1
   uml2 = im1file->BY473EM_DIFF_ATNJJA0_V_UASD0JF1
   vml2 = im1file->BY473EM_DIFF_ATNJJA0_V_VASD0JF1
   smr2 = im2file->BY473EM_DIFF_NINO34JJA0_V_TSJJA0
   umr2 = im2file->BY473EM_DIFF_NINO34JJA0_V_UASJJA0
   vmr2 = im2file->BY473EM_DIFF_NINO34JJA0_V_VASJJA0
   sml2 = sml2 + lndmsk
   uml2 = uml2 + lndmsk
   vml2 = vml2 + lndmsk
   smr2 = smr2 + lndmsk
   umr2 = umr2 + lndmsk
   vmr2 = vmr2 + lndmsk
;=== third data ===;
   i3file = addfile(idir(1)+fileNames(4),"r")
   i4file = addfile(idir(2)+fileNames(6),"r")
   i5file = addfile(idir(2)+fileNames(8),"r")
   i6file = addfile(idir(1)+fileNames(5),"r")
   i7file = addfile(idir(2)+fileNames(7),"r")
   i8file = addfile(idir(2)+fileNames(9),"r")
   sl3 = i3file->RMAP_ATNJJA0_V_TSD0JF1
   ul3 = i4file->PVAL_RMAP_ATNJJA0_V_UASD0JF1
   vl3 = i5file->PVAL_RMAP_ATNJJA0_V_VASD0JF1
   sr3 = i6file->RMAP_NINO34JJA0_V_TSJJA0
   ur3 = i7file->PVAL_RMAP_NINO34JJA0_V_UASJJA0
   vr3 = i8file->PVAL_RMAP_NINO34JJA0_V_VASJJA0
;--- observation mask ---;
   sl1 =  where( abs(sml1).ge.0, (sl1 + lndmsk), lndmsk)
   sr1 =  where( abs(smr1).ge.0, (sr1 + lndmsk), lndmsk)
   sl2 =  where( abs(sml2).ge.0, (sl2 + lndmsk), lndmsk)
   sr2 =  where( abs(smr2).ge.0, (sr2 + lndmsk), lndmsk)

   ul1 =  where( abs(uml1).ge.4, (ul1 + lndmsk), lndmsk)
   ur1 =  where( abs(umr1).ge.4, (ur1 + lndmsk), lndmsk)
   ul2 =  where( abs(uml2).ge.4, (ul2 + lndmsk), lndmsk)
   ur2 =  where( abs(umr2).ge.4, (ur2 + lndmsk), lndmsk)

   vl1 =  where( abs(vml1).ge.4, (vl1 + lndmsk), lndmsk)
   vr1 =  where( abs(vmr2).ge.4, (vr1 + lndmsk), lndmsk)
   vl2 =  where( abs(vml1).ge.4, (vl2 + lndmsk), lndmsk)
   vr2 =  where( abs(vmr2).ge.4, (vr2 + lndmsk), lndmsk)

   ul3 = (ul3 + lndmsk)
   vl3 = (vl3 + lndmsk)
   ur3 = (ur3 + lndmsk)
   vr3 = (vr3 + lndmsk)

delete([/sl1@long_name,sr1@long_name,sl2@long_name,sr2@long_name,sl3@long_name,sr3@long_name/])
;--- output variables ---
  diro    = "./"
  pltName = "fig03"
  outFile = diro+pltName+".nc"
  system("\rm "+outFile)
  fout         = addfile( outFile,"c")
  fout->ul1 = ul1
  fout->vl1 = vl1
  fout->sl1 = sl1
  fout->ul2 = ul2
  fout->vl2 = vl2
  fout->sl2 = sl2
  fout->ul3 = ul3
  fout->vl3 = vl3
  fout->sl3 = sl3

  fout->ur1 = ur1
  fout->vr1 = vr1
  fout->sr1 = sr1
  fout->ur2 = ur2
  fout->vr2 = vr2
  fout->sr2 = sr2
  fout->ur3 = ur3
  fout->vr3 = vr3
  fout->sr3 = sr3



;--- plotintg ---
   wks  = gsn_open_wks("png","fig03_newcolor") ; send graphics to PNG file
opt = True
opt@Slat = -40
opt@Nlat =  40
opt@Elon = 100
opt@Wlon = 440

opt@expName = ""
opt@verName = ""
   plot = new(6,graphic)
   plot(1) = plot_shade_con_vec_2Dpattern(wks,ul3,vl3,sl3,tsmk1,opt) ; HadISST
   plot(3) = plot_shade_con_vec_2Dpattern(wks,ul1,vl1,sl1,sml1,opt) 
   plot(5) = plot_shade_con_vec_2Dpattern(wks,ul2,vl2,sl2,sml2,opt) 

   plot(0) = plot_shade_con_vec_2Dpattern(wks,ur3,vr3,sr3,tsmk2,opt) ; HadISST
   plot(2) = plot_shade_con_vec_2Dpattern(wks,ur1,vr1,sr1,smr1,opt) 
   plot(4) = plot_shade_con_vec_2Dpattern(wks,ur2,vr2,sr2,smr2,opt) 

;---Draw and panel all four plots
MString = "        Observed and simulated interactions between the    ~C~   tropical Pacific Ocean and the equatorial Atlantic Ocean ~C~ ~Z55"
  pres             = True
  pres@gsnMaximize = False
  pres@gsnPanelLabelBar   = True
  pres@pmLabelBarWidthF   = 0.8
  pres@lbLabelFontHeightF = 0.015
  pres@gsnPanelMainString        = MString
  pres@gsnPanelMainFontHeightF   = 0.016
  pres@gsnPanelFigureStringsFontHeightF = 0.012
  pres@gsnPanelFigureStrings     = (/"a)","d)", "b)","e)", "c)","f)"/)     ; add strings to panel
  pres@gsnPanelFigureStringsJust = "TopLeft"                      ; new resource added in NCL V6.4.0
  gsn_panel(wks,plot,(/3,2/),pres)

end
